#对颗粒相进行轴向处理
#by zyp
import numpy as np
import pandas as pd
import matplotlib
import matplotlib.pyplot as plt
import math
print('cal...')
n = 19
name_data = 'par34_33_300s'

name = np.array(pd.read_csv(name_data, nrows=n + 1))
# 读取dat文件数据
real_data = 2  # 从第j个开始是物理数据不要改
look_data = 7  # 想看哪个数据7是c，11是times
species = 1  # 待生剂编号
data = np.loadtxt(name_data, skiprows=n + 2)
total_rows = data.shape[0]
total_cols = data.shape[1]
data = data[data[:, n - 2].argsort()]

####轴向######
y_min = min(data[:, n - 2])
y_max = max(data[:, n - 2])

point_number = 50  # 高度分几份

interval = (y_max - y_min) / (point_number - 1)
result = np.zeros((point_number, n))
result[:, 0] = np.linspace(y_min, y_max, point_number)
interval_use = interval / 2

j = real_data - 1
while j < n - 3:
    print(name[j + 1])
    i = 0
    k = 0
    sum = 0
    count = 0
    while i < total_rows:
        if data[i,4]!=species:
            i+=1
            continue

        if result[k, 0] - interval_use < data[i, n - 2] <= result[k, 0] + interval_use:
            sum = sum + data[i, j]
            count += 1
            i += 1
        else:
            if count == 0:
                count = 1
                print('count=0')
            aver = sum / count
            result[k, j] = aver
            k += 1
            sum = 0
            count = 0
    result[k, j] = aver
    j += 1

np.savetxt(name_data + '.csv', result, delimiter=',')
plt.plot(result[:, look_data - 1], result[:, 0])
plt.title(name[look_data])
plt.show()
print('finish')